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semileptonic electroweak interactions 
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Scmi-leptonic electroweak interactions in nuclei - such as j3 decay, n capture, charged- and neutral- 
current neutrino reactions, and electron scattering - are described by a set of multipole operators 
carrying definite parity and angular momentum, obtained by projection from the underlying nuclear 
charge and three-current operators. If these nuclear operators are approximated by their one-body 
forms and expanded in the nucleon velocity through order \p\/M, where p and M are the nucleon 
momentum and mass, a set of seven multipole operators is obtained. Nuclear structure calculations 
are often performed in a basis of Slater determinants formed from harmonic oscillator orbitals, 
a choice that allows translational invariance to be preserved. Harmonic-oscillator single-particle 
matrix elements of the multipole operators can be evaluated analytically and expressed in terms of 
finite polynomials in q 2 , where q is the magnitude of the three-momentum transfer. While results for 
such matrix elements are available in tabular form, with certain restriction on quantum numbers, the 
task of determining the analytic form of a response function can still be quite tedious, requiring the 
folding of the tabulated matrix elements with the nuclear density matrix, and subsequent algebra 
to evaluate products of operators. Here we provide a Mathematica script for generating these 
matrix elements, which will allow users to carry out all such calculations by symbolic manipulation. 
This will eliminate the errors that may accompany hand calculations and speed the calculation of 
electroweak nuclear cross sections and rates. We illustrate the use of the new script by calculating 
the cross sections for charged- and neutral-current neutrino scattering in 12 C. 



I. INTRODUCTION 

A common task in nuclear physics is to connect an observable - a rate for decay, the diffraction pattern seen in 
inelastic electron scattering, etc. - back to the underlying nuclear structure physics. While this can always be done 
implicitly through numerical calculations, one of the attractive properties of the harmonic-oscillator shell model is 
that this connection can be made analytically, at least in the case of one-body observables That is, each observable 
can be expressed as a simple function involving finite polynomials in q 2 , the square of the three- momentum transfer to 
the nucleus, with the polynomial coefficients depending on the one-body density matrix. This allows one to determine 
quickly how a given experimental result constrains the underlying nuclear physics (the one-body density matrix), 
and to dentify the unconstrained degrees of freedom in the density matrix that will enter into predictions of other 
observables, and thus will have to be taken from models. 

Specifically, electroweak processes such as j3 decay, muon capture, charged- and neutral-current neutrino reactions, 
electron scattering, and photo-absorption and 7 decay are traditionally described as one-body processes, possibly with 
corrections added to account for exchange currents, final- or initial-state distortions of the lepton wave function, etc. 
That is, the charge and three-current nuclear operators are expressed as the sum of the currents of the constituent 
nucleons. As described below, the many-body matrix elements of such operators can then be expressed as sums over 
single-nucleon matrix elements. Some time ago harmonic oscillator matrix elements - that is, the coefficients of the 
polynomials mentioned above - were tabulated by Donnelly and Haxton |2|, though with some restrictions to reduce 
the size of the compilation. The tables include single-nucleon states up to the closed shell at proton/neutron number 
126, e.g., all oscillators shells for principal quantum number N = to 5, as well as the j = 13/2 subshell for N=6. 
(N determines the harmonic oscillator energy, (N + 3/2)frui.) The tables also restricted transitions to \Nf — iV,| < 2. 

In today's context, the main shortcoming of the tables is not these restrictions, but the need to transcribe the 
results when using them in calculations. This is tedious and provides opportunities for errors. Thus here we provide 
a Mathematica 3] script to generate equivalent results, without restrictions on the choice of quantum numbers. 
This script will allow users to efficiently manipulate results within Mathematica, so that cross sections, rates, and 
relationships between various observables can be obtained quickly and without error, in closed form. 
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II. ONE-BODY OPERATORS AND THE DENSITY MATRICES 



We consider electroweak nuclear multipole operators Oj-t that are one-body, e.g., have the first-quantized form 
Oj ; t(*), and carry definite parity and angular momentum. This operators, which are formed from the single- 
nucleon currents and take into account the spatial extent of these currents, can then be evaluated in terms of the 
one-body density matrix, 



■l:: Tj[\() J:r \ll,: T, £ ^£(|a|, |/3|)<|a|; #-.0 J{r : \\p\- t ±). (1) 

Here (J;T), (Jf,Tf) and (J^Ti) are the angular momentum and isospin of the operator and of the final and initial 
nuclear states, respectively, |a| = {n a (l a ^)j a } represents the set of nonmagnetic single-nucleon spatial and spin 

quantum numbers, and :: indicates a matrix element reduced in both angular momentum and isospin, i.e., 
( J,M,; T,M T , \0j. :T \ J,Mi : T,M T , ) = 

i-n"-"'^-"" ( _% f I, * ) ( X, Mr M T , ) (J,:T,id J;T ij,^. (2) 

The sums extend over complete sets of quantum numbers |a| and |/3| - which would normally correspond to the 
single-nucleon basis used in constructing the Slater determinants for the initial and final nuclear states. 

For a one-body operator O, Eq. ^ is exact and factors the operator physics - embodied in the single-particle 

matrix elements (\a\; h'.'.O j ; t"-\/3\; \) - from the many-body nuclear physics, which is contained within the one-body 
density matrix 

1>jtW) = [jpj< J / ;T / :: [ c H ® 5 I/3|Ut!^;T 4 ). (3) 

Here the single-nucleon creation and annihilation operators carrying good angular momentum and isospin are defined 
by 



|a|mo, ;m,{ 

{ _ 1)ja - ma+ X/2- mtaCiai _ m ^_ m ^ (4) 



where m a and mt a are the single-nucleon magnetic quantum numbers for angular momentum and isospin and a — 
{\a\m a ;m ta }- 

In practice, of course, the density matrix of Eq. ([3]) is often not calculated from first principles: for the classi- 
cal nuclear physics problem of bound nucleons interacting via a two-body (or two- plus three-body) potential, the 
Schrordinger-equation many-body problem may not be tractable. Instead, wave functions are taken from nuclear 
models, such as the shell model, where model assumptions restrict the single-particle states that can be occupied. 
An effective interaction, often determined empirically, is employed to account for the effects of the omitted degrees 
of freedom. This leads to a finite basis of Slater determinants for the many-body problem and thus a Hamiltonian 
problem that can be solved by direct diagonalization. 

To take the example we will use in this paper, consider the case of a Ohui shell-model calculation for 12 C. That is, the 
included space consists of all Slater determinants where the Is shell is fully occupied and eight valence nucleons are 
distributed in all allowed ways among the twelve single-particle orbitals within the lp shell. Initial and final nuclear 
states determined from such a calculation would both carry positive parity. Under this assumption the density matrix 
for a transition between two such states simplifies, 



{lpjT(a, /?)} -> {'0jt( 1 Pi/2: lPl/2), "0JT + ( 1 P3/2, IPl/l), VjT (^3/2) lPl/2)) V'jt( 1 P3/2, I-P3/2), ^jt( 1s 1/2, lSl/ 2 )}; 
WjT (!P3/2) 1 Pl/2) = ^7^(^3/2, lPl/2) ± - 0jt( 1 Pi/2> ^3/2); ^jt( 1s 1/2, lSl/2) = 25f^S Jfi 5 T ,Q (5) 

where ipjT ± (^P3/2, 1.P1/2) ar e defined in anticipation of the time-reversal properties of the operators we will discuss 
below. Note that the lsi/ 2 shell only contributes to J=0 and T=0 elastic transitions. 

So far the single-particle basis has not been specified: in many cases it can remain unspecified through the calculation 
of the one-body density matrix. This would be the case in shell-model calculations where matrix elements of the 
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effective interaction are treated as parameters, determined by fitting energy levels, e.g., as Cohen and Kurath [4( did 
in the lp shell. Such a shell model calculation would yield specific numerical values for the density matrix elements 
of Eq. (HJ). However, when Eq. (pj is invoked to evaluate matrix elements, generally a single-particle basis must be 
specified. 

The harmonic oscillator basis is an attractive choice for this purpose because the single-particle matrix elements for 
electroweak interactions can be evaluate analytically and have a simple form involving polynomials in the magnitude 
of the three-momentum transfer to the nucleus. This is very useful for determining the functional form of electroweak 
response functions. Even if the density matrix is unknown, response functions can be expressed in a form where the 
density matrix elements are parameters. This allows one to impose experimental constraints on model calculations - 
to quickly determine what degrees of freedom in the density matrix are already constrained by experiment, so that 
these can be eliminated, reducing the model dependence in predicting other process. In an example we will discuss 
later, charged- and neutral-current neutrino reactions in 12 C involving the triad of JT = 1+1, Mt = (—1,0,1) 
excited states, the density matrix elements are sharply constrained by results from electron scattering (and gamma 
decay), (3 decay, and /i capture In some truncated model spaces there may be a sufficient number of experimental 
constraints to determine all of the contributing density matrix elements: this would then free one from dealing with 
nuclear models. 



III. THE BASIC NUCLEAR OPERATORS AND THEIR HARMONIC OSCILLATOR MATRIX 

ELEMENTS 

The semileptonic weak nuclear operators Ojm-,tm t °f interest in this paper are one-body and can be expressed as 
a product of space-spin and isospin operators. In first quantization 



A 

E 



Oj M (i)lTM T (i), 



(6) 



where the space-spin operator Oj includes an isospin label because it contains couplings, like magnetic moments, that 
have different strengths depending on the isospin, e.g., the isoscalar magnetic moments is not equal to the isovector 
magnetic moment. Semileptonic nuclear cross sections and rates can be expressed in terms of nuclear matrix elements 
reduced in angular momentum 



(Jf\\dj,TM T U 

(-1) 7 ' - 



(J f ; T f M Tf 1 1 Oj,tm t \\Jr, T z M Ti 
ip f j'T(\a\, 



r, \i r , I Tf T Ti \ 



-M Tj M T M Ti 



\Oi 



9.\\It\\o) 



(7) 



where the arrow in the first line indicates we are ignoring small Coulomb and other charge-independence-breaking 
terms, so that nuclear states can be labeled as eigenstates of good isospin (a common assumption in nuclear structure 
calculations). The second line follows from Eq. ([I}. 

The one-body space-spin operators Oj governing semileptonic electroweak interactions are derived in standard 
references [7J by expanding the vector and axial- vector charge and current operators through order \p\/M. They can 
be expressed in terms of the seven single-particle operators 



Af'(gx) 



A'f'(gx) 



Mfy (<pc) • -v 



-Vx M 



M / 
J.J 



-J 1 ' 2 Mfj +1 (gx) + (J + l) 1 / 2 M¥jU (gx) 



Mj 
J J 



E'/^gx) ee -i 



ivx Mf/ (gx) 
g 



O^(gx) 



'(gx) 



-VMf J (gx) 
g 



a = [J]- 1 M# +1 (gx) + ( J + l) 1 / 2 MfjU (gx) 

3 = [J]- 1 \(J + I) 1 ' 2 Mff +1 (gx) + J 1 ' 2 M^_ x (gx) 



Mf J (gx) a - -V, 

q 



(8) 
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where q is the magnitude of the three-momentum transferred to the nucleus and [J] = \J2J + 1. Because of operator 
time-reversal properties it is helpful to replace O by a new operator 

n^'(gx) = n^'fox) + l^ M '(qx). (9) 

The multipole operators in Eqs. (J5J) are constructed from spherical Bessel functions, spherical harmonics, and vector 
spherical harmonics Q 

Mf'(gx) = JJ (qx)Y J M -'(n x ) 

M% (gx) = j L (<7*) (Q x ) (10) 

where 

(O x ) = [Y L ® ei } JMj = 5Z<£mlA| (Ll)JMj)Y^(Q x ) e 1A . (11) 

The spherical unit vectors eiA are defined by ei±i = T{ e x ±i e j/)/v / 2 and eio=e z . 

The single-particle basis |a) = |fta(^aT/2)ja m a) of Eq. (TT|) has the general coordinate-space form 

^UaW^^l^i/jk™. ( 12 ) 

where £]y 2 is the Pauli spinor. The angular and spin portions of the single-particle operator matrix elements can then 
be performed, leaving only matrix elements between unspecified radial wave functions, 

{n\l>\)j>\\Mj{ q *)\\n{l\)j) = -!=(-) p/j {l] [f] {j] [j] { ^ J / 3 } ( ) {n ' 1 ' 3 ' ^WW) 

{n\l'l)f\\ M, 7L (qx) ■ a\\n(ll)j) = {^(-f [I] [j] [L] W j | I 1 } ( q q ) {n ' l ' j ' \^ nl ^ 



V A' 1 



n'{l'\)j'\\ M, 7L fox) • iv|KJ±)j) = J=(-)WV2 p/j y/j [£] [j] j ? J j 



x {-(/ + D 1 / 2 [I + 1] { f J , ; 7 ! } ( I L Q 1 + 1 ) W |iz(p) (| - |ny) 

(n'(Z'|)/||M.,(gx)(? • iv||n(J±)j) = ^=(-)'' [0 [?'] [?] [2j - i] [J] { J 2 /l | } } ( q q ^o" ' ) 

x|-^ m/2 (n'6-' |j./(p) - + S jtl _ 1/2 {n'l'j' \jj{p) (± + |nJj) J. (13) 

Here the dimensionless coordinate p = ga;, while the radial matrix elements are defined by 

(n'l'j'\9(p)\nlj) = [ x 2 dxR* n , Vl ,(x)9( P )R nl] (x) (14) 



for 



jj(p) 



8( P ) = I Jj(p) (±- L P 

Up) ( 



(15) 



_d_ , |+1 
dp "r" p 



To complete the evaluation a specific choice for the single-particle radial wave functions must be made. The 
harmonic oscillator is an attractive choice, convenient both for nuclear structure structure reasons (for example, in 
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certain model spaces spurious center-of-mass degrees of freedom can be removed exactly) and because the radial 
matrix elements of Eqs. (|15[) can be evaluated analytically. The properly normalized radial wave functions are 



Rnlj(x) = R n l{x) 



2e z 



b 3 (n~ l)W(n + l + ±y+! 



- 1/2 d n-l 






' z n+l-l/2 e -z 


dz 11 - 1 





(16) 



where z = (x/b) 2 and b is the oscillator size parameter. We employ a notation where the nodal quantum number 
n = (N — l)/2 + 1, where the principal quantum number N = 0,1,2,... Thus n= 1,2,3... A given harmonic oscillator 
state can be labeled by n,l,j or, equivalently, N,j: as N determines the parity of the shell, knowledge of N and j 
uniquely determines I. Our Mathematica script uses N,j to identify single- nucleon states. 

The radial integrals appearing in Eqs. (| 13[) can be evaluated analytically for harmonic oscillator states, 



(n'l'\ 3L (p)\nl) = (2 L + m V L/ *e- V j(P W - l)!T(n' + I' + |)r(n + 1 + ^^ 
^ '" fe=0 fe'=0 

1 r[i(/ + /' + L + 2fc + 2fc' + 3)l rl/ , , , 

(n'iWxJ; - = (2 I + i)!! y(£ " 1)/2e " V (n ~ 1)!(n ' ~ 1)!r(n ' + + ^ )T{n + l + ^ 
^^ (_i)fc+fc' 1 r[i(x + / + r + 2fc + 2fc' + 2)] / 

* hvho klk ' 1 (n-l-k)\(n'-l-k')\ T[l + k+l}r[l' + k' + i} \ 1 ( + + + + + ' 
x i-Fi[i(i - Z - I' - 2k - 2k' - 1);L + + (2k) i-Fi[|(i - I - I' - 2k - 2k' + I); L + §; y] j (18) 

<^U(p)(|; + ^)M> = + i)!! y(I " 1)/2e " - 1)!(n ' - 1)!r(n ' + + ^ )r(n + ' + 



(2£ + l)! 

^ ^ (_i)fc+fc' 1 r[^(x + i + r + 2fc + 2fc' + 2)] j , . 

X hiikn km {n-l-k)\{n'-l-k')\ T[l + k + |]r[/' + k' + f] 1 ^ +L + 2k + 2k +2) 



x i-Fi[|(i - i - Z' - 2fc - 2Jfe' - l);i + + (2/ + 2fc + 1) - Z - l' - 2k - 2k' + 1);L + §; y] j (19) 

where y = (qb/2) 2 . Here 1F1 is the confluent hypergeometric function 

In the present application a is a nonpositive integer, so that this series is a polynomial of order —a in y. 

Thus the single-nucleon harmonic-oscillator matrix elements of the seven basic operators are completely determined 
in terms of simple functions. One finds 

(n'{l'±)j'\\Tj(q*)\\n{l\)j) = -Ly^)/^-^) , (21) 



where K = 2 for the normal parity operators M, A', S; K = 1 for the abnormal parity operators A, £', f2 (or O'). 
The tabulations of Ref. [2| provide the polynomials p(y) for such matrix elements, with the restrictions described 
previously. The Mathematica script accompanying this paper generates the full expression in Eq. (|2ip . without 
restrictions on the quantum numbers and in a form that can be easily manipulated within Mathematica to generate 
analytical expressions for form factors, rates, etc. 

An important property related to the time-reversal properties of these operators is 

(n>(l'±)j>\\Tj( q x)\\n(li)j) = (-l)\n(l^j\\T,( q ^\\n'm)j') (22) 

where A = j' — j for the operators M, A, £', and £", and A = j' + j for the operators A', S, and fi'. The 
operator £1' of Eq. ^ was defined because f2 lacks such a simple transformation property. This "turn around" 
property is important because it leads to a dependence of observables of density matrix element combinations such 
as ipfy (lp 3 /2,lp 1/2 ). 
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IV. RELATIONSHIP WITH SEMILEPTONIC PROCESSES IN NUCLEI 



To illustrate the use of the Mathcmatica script, we will consider various charged-current and neutral-current neutrino 
reactions from the J^T = + ground state of 12 C: 

12 C(g.s.)( Vei e~) 12 N(g.s.) 12 C(g. S .)(D e ,e+) 12 B(g. S .) 
12 C(g.s.)(v, i/) 12 C(15.11 MeV) 12 C(g.s.)(P, ^) 12 C(15.11 MeV) 

12 C(g.s.)(is, v) 12 C(\2.1 MeV) 12 C(g.s.)(D, v) 12 C{12.7 MeV) 

The first four transitions are isovector, populating the triad of J + T = 1 + 1 states formed by the ground states of 12 N 
(Mr = 1) and 12 B (M T = -1) and the 15.11 MeV (M T = 0) isovector analog state in 12 C. The last two are isoscalar, 
populating the J + T = 1+0 12.7 MeV state in 12 C. The four-momentum transfer for these reactions is defined by 

<f = P?-Pjf = k»-ki; = {w,(fi (23) 

where P-* , P?, k£ 7 and are the four-momenta of the initial nucleus, final nucleus, incident neutrino, and scattered 
lepton (electron/positron or neutrino/antincutrino). The charged-current cross section in the extreme relativistic 
limit can be written ll| 



fda_ 



oo 



-G 2 F cos 2 9 C F(Z, e) —— cos 2 I ( Jf\ \M j + \ J, 

ir Ui + 1 2 '-J q 
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=F 2 tan - 



tan 2 — 



2 i |/7| \ /f-mag i 



J=l 
.2 



-4+ tan - 
q A 2 



E[K J /ii^ii^>r + i(^/ii^ m 

Y,Re((j f \\Tj ma9 u)(jf\\v l u: 



,7=1 



(24) 



where Gf cos 9 C is the weak coupling constant, 9 the angle between the incident neutrino and outgoing electron, 
q 2 = uj 2 — q 2 , and e is the energy of the outgoing electron/positron. The function F(Z, e) corrects the electron phase 
space for the effects of Coulomb distortion in the field of the daughter nucleus of charge Z. This expression assumes 
e >> m e and neglects nuclear recoil. 

The operators Mj, Cj, T el , and T mag are the familiar charge, longitudinal, transverse electric, and transverse 
magnetic projections of the weak charge-changing hadronic current, 



(25) 



where JTJ, i = 1,2,3 are the three components of an isovector. The weak current is made up of vector and axial- vector 
components 



qt _ n 



J, 



(26) 



Multipole projections of this current are taken to exploit the angular momentum and parity quantum labels that 
are normally carried by the initial and final nuclear states. The multipole operators are defined by 

Mj M j;TM t = [dxMf J (qx)jo(x) T M T =Mj M j;TM T +Mj Mj . TMT 
£jMj;TM T 



2 JMj;TM T 



1 JMj-TM T 



q 

- I <Ix 

q 



VMj /j (qx) • J(-x)tm t - Ljaij-tMt + L 5 JMj . TMt 
Vx Mfj (gx)] • J( X ) TAlT = Tf Mj . TMT +Tf! Ij . TMT 
dx Mfj (gx) -J{-x)tm t - T7m.,:Tm t + ^?m%m t 



(27) 



where we have separated the vector and axial- vector (with superscript '5') contributions in Eqs. (|27[) . The operators 
Mj, Lj, ff, and f] la95 are normal parity operators, that is, Avr = (-1)' 7 . The operators M% L 5 j, ff 5 , and fj ag 
are abnormal parity, Air — (— 1) J+1 . Also note that, for a conserved vector current, Lj can be eliminated as 

to 



L 



JMj-TMi 



-M 



q 



JMj-TMt- 



(28) 
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Although contributions from two-body (exchange) currents can be included, more commonly the nuclear charges 
and currents are modeled as the sum of the one-body contributions from the constituent protons and neutrons. For a 
single free nucleon one has the following general forms for the matrix elements of the vector and axial-vector isovector 
currents for charge-changing weak interactions 



(fc'V; ±m r |J±(0)|feA; \m t ) = u(k> '\'){f[ 1] '(<£) 7m - iF^\ql)a llv q v }u{k\){\m v \r ± \\m t ) 
(fc'A';imHJf (0)|fcA;im f ) = u(k' \'){F^ (ql) W - F^\ql) l5 q^u{k\){\m^r ± \\m t ) (29) 

where we have omitted second-class scalar (which violates conservation of the vector current) and axial-tensor cou- 
plings. [Our notation in this paper follows Bjorken and Drell [7|], including gamma matrices, metric, and current 
definitions, with the exception of the definition of the form factors, which follows that of Ref. 

A standard non relativistic reduction of the matrix elements through order 1/Mjv, in which momenta are interpreted 
as gradients operating within the nucleus, leads to expressions for the multipole operators of Eqs. (f2"7| in terms of 
the seven basic single-particle operators defined in Eqs. ©: 

Mf Mj (qx) = F[ l \ql)M™J{qx)T ± 



Tf M %*) = ^(F?\ql)A'f'(qx) + ^( q l)^( qX ))T ± 
T-^x) = -^ { FP(ql)A^(q X )-^\ql)^(q X ))r ± 
Mj±>x) = ^(^(^ ) ^( ?x) + | w F«(g2)E^( g x))r ± 



N 



L 5 ^M) = l (id 1) (^)-^^ ) (^))s:;^(gx)r ± 



T™f ^x) = F£\<t>)J%'(qx)T ± (30) 



'■ JM. 

where M N is the nucleon mass and /i (1) (<?^) = F [ 1] {q 2 ^) + 2M N F^\ql). 

The nucleon isovector form factors appearing in Eq. ([3D]) include the charge form factor (.F^ (0) = 1), the magnetic 
moment (/i^^(0) ~ 4.706), and the axial form factor (F^ (0) ~ —1.26). Various dipole and other representations of 
the momentum dependence of F± , F^\ and F^ are available; the momentum dependence of the pseudoscalar form 
factor Fp^ is found to be 

4V„)~^P, (3D 

assuming pion-pole dominance and the Goldberger-Treiman relationship The isospin raising/lowering operator 
t± is related to the spherical projections of isospin Pauli operator r\ m by 

4 =1 = ^Ut ± =T^t 1±1 . (32) 



Charged- current neutrino scattering off 12 C: Now we consider the example of charged-current neutrino scattering off 
the JT = + 12 C ground state, leading to the 1+1 grounds states of 12 B and 12 N. We treat these states as isospin 
eigenstates, ignoring isospin breaking due to electromagnetic or other charge-dependent interactions. Thus Eq. (pM)) 
simplifies because of the angular momentum and parity constraints imposed by the + — > 1 + 1 transition and because 
the matrix elements can be reduced in isospin, 
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3tt 
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tan 
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where a factor of 1/3 results from isospin reduction, 



(33) 



(l;l±l||O 1;1±1 ||0;00) = -= 1: 1 \'.0\- \ ::(): (I . 



(34) 



Eqs. (j3U|) and (|3ll|) thus determine the cross section, once the doubly-reduced single-particle matrix elements are 
evaluated. Eqs. (JJ) and (J32J) yield 



(i+;iiA il ;;o+ ; o)= Yl 4ii-,T=i(HW\)(H,^l\ 



H,\0\ 



£ ^i(l«|,|i9|)<H IIO^II |/3|)(T^)(i||r||i)=TV3 ]T ^(\a\, |/3|)(|a| ||0 

l«|,|/3| l«|,|/3| 



(35) 



where we have used (|||r|||> = V6. The remaining spin-spatial matrix element can then be evaluated with the 
Mathematica program, by virtue of Eq. (|30|) . As discussed early, we truncate the density matrix to the lp shell 
and assume a harmonic oscillator basis, but otherwise keep the elements of the density matrix arbitrary: In the next 
section we will describe how the Mathematica code can be used to provide the needed matrix element, once this 
truncation has been made. The result is an analytic expression for the cross section that encompasses any shell-model 
calculation restricted in this way, 



™C(u B ,e-rN(g. s .) = ^2°f cos2 WZ, e) cos 2 °-e 2y 
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F, 
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2 /, „ \ , 4 V2, y N , 2\/l0, 4 . , 
--(1 - 2y)1> 11 - -^-(1 - |)V-3i + -^-(1 ~ 5^33 



=F2 tan • 



e 



1 1/2 



2 1 tan ' o 
7 2 



d-|)^ + ^(l-^33 U 



,-.(!) 



2 \/2 _ "\/l0 

3^11 + ^31 + ~ 



33 



^H(l-Wn-^(l-|WS + ^(l-^33 



(36) 



where it is understood that all single-nuclcon couplings are described by form factors, e.g., F^ — F^\qf t ). In this 

expression we used the short-hand notation for the density matrix of tf>2j a 2j p = 'i'jZX t=i° °(l a l> e ^ c - That is, 
the final and initial nuclear state labels /, i are suppressed, as are the J, T of the transition, since these quantum 
numbers are uniquely determined for the transition ( J = 1, T = 1). 
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Neutral current scattering to the 1 + 1 15.11 MeV state of 12 C: The hadronic neutral current is 

J"° = J*-28m 2 6 w J™ (37) 

where J 3 is third isospin component of the weak current of Eq. (125[) . J^ m is the electromagnetic current, and Qw the 

Weinberg angle. Because of the presence of , the neutral current includes an isoscalar contribution. Eqs. (f!?5|) 
and (|37[) allow one to compare the isospin dependence of charged and neutral-current isovector transitions, 

7y =1 = t± <-> It=i = [1 — 2 sin 2 8 W ]^, vector current 

«-> It=i 5 — axial — vector current (38) 
The reduction in isospin between nuclear states (compare Eq. (|34|) ) yields 

(1; 10||O 1;10 ||0; 00) = -^=(1; l:!6i.i::0; 0), (39) 
while the single-particle matrix element for T3/2 (compare Eq. ([55)) ) yields 

~ V^(M,l/?l)(M WO^W \(3\){3\\tU) = ^§ £ ti-A(n\0\)(\<x\ \\0P\\ |/3|) • (40) 
\*\M M,l/?l 

Thus there is an overall factor of 1/2 in the cross section due to isospin, relative to charged-current result in 12 C. 

Finally, the cross section for neutral-current scattering of neutrinos off nuclei can be obtained from Eq. ([24| by 
replacing e with e„, the energy of the scattered neutrino, and by removing F(Z, e) and cos 2 9 C . It follows 



da " 



dfl ) „c(i+i) 4 ^ 2 



^Gi cos 2 --e~ 2y 



2 c(p,p') 12 c(i + i) 



M A q [ ' 



(i) _ 1 



■^(1 + 2y)i> n - i(l - y)^ + ^(1 - ^33 



2?2 +tan 



V2 



tf> ( -|d - - ^(1 - + ^(1 - ^33 



+ ^(l-2sin 2 ^) 2 x 



if > ( hn + ^sl + ^33 ) + ^ ( ^ 2^ - ^(1 - §)fe + - 



=F2tan 



3 r±± 3 

? 
9 



%+tan 2 ^ 



3 

1/2 



4^(1 - 2sin 2 6 W ) (-1(1 - 2y)^ - ^(1 - f)^ + ^(1 - | s )fc, ) x 



V2 VTo 



,(!) 



4v^,, ,_ 2\/To. 



33 



(41) 



where, of course, both the Coulomb correction and effects of the Cabibbo angle have been removed. 

Neutral current scattering to the l + 12.7 MeV state 12 C: The transition to the l + 12.7 MeV state is generated by 
the isoscalar neutral current which, according to Eq. (|37[) . involves 



It£q = — (2 sin 2 W )—, vector current 
Ix=q 5 = 0, axial — vector current 



(42) 
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That is, the isoscalar neutral current is purely vector which, for example in 12 C, thus limiting the contributing 
multipoles to Tj 1 ^ . The isospin factors are trivial 



(l;00||O 1;00 ||0;00) = (1; 0::Oi i0 ::0; 0), 
while the single-particle matrix element for the unit isospin operator is 



\<*\,\P\ 



\<*\M 

It follows that the cross section is 



H,\0\ 



(43) 



a\ NO 



(i) 



(44) 



da 

p(0) 



12 C(u,u') 12 C (1+0) 
12 C(P,P') 12 C(l + 0) 



7^1 



^31 



G'i sin 4 W cos 2 -e"^ 



2g 2+tan 



q 

M 2 



10 



.S3 



/' 



--(1 - 2y)xhi - -^-(1 - |)^3i + -^-(1 - 



1 \ 3 r 3 3 r / 2 

The isoscalar vector couplings are functions of g 2 , with (0) = 1 and //(°) (0) = 



(45) 



V. REDUCED MATRIX ELEMENT FUNCTIONS IN MATHEMATICA 

The Mathematica package SevenOperators has been developed and made available online Q to replace and extend 
the tabulated results of Ref. Q. The main motivation for this package is to make calculations of semileptonic 
electroweak cross sections and decay rates far less tedious: users can couple it to cross section and rate formulas, 
explore various truncations of the density matrix, determine density-matrix constraints imposed by known rates or 
cross sections, etc. Such closed-form expressions previously would have required tedious hand calculations. The 
SevenOperators package has been cross-checked against Ref. 0, including the production and posting Q of tables in 
the form of Ref. . 

In this section we give a basic description of SevenOperators, and provide useful technical notes. As an example of 
usage, we also show its application to the calculation of the 12 C transitions, which were illustrated in the previous 
sections. 

The package is available in two forms. The first is a notebook version (7operators.nb) that was developed on 
Mathematica 5.2. The second is a modular version (package. tar. gz) suitable for use with both the graphical and 
text-based versions of Mathematica. It consists in a .tar.gz archive containing a master file (7o_master_file.m) and 
several input files, with self-explanatory names, each of them properly commented. The whole package is run by 
simply running the master file. 

SevenOperators returns seven functions, which give the reduced matrix elements of the seven basic operators between 
single-particle states, in the form of Eq. pTjl . The functions have the variables y, N' , f , N, j, J as arguments, and 
their Mathematica syntax is given in Tables [J and [TT] for normal and abnormal parity operators respectively. 



TABLE I: The functions describing the reduced matrix elements of the normal parity operators, as implemented in SevenOp- 
erators Mathematica package For reference, next to each function we give the reduced matrix element evaluated in Seven 
Operators (see Eqs. (01). Note that the labels N and j fully specify the standard space-spin nonmagnetic quantum numbers 
for a harmonic oscillator, n(l l/2)j. 



Mathematica function 


Reduced matrix element 


Mj[y,{N', j'},{N, j}, J] 


(n'(l' l/2)j||AC(gx)||n(Z l/2)j) 


SigmaJ[y, {n', j'}, {N, j}, j] 


<n'(J' l/2)f\\Z^(qx)\\n(l 1/2) j) 


DeltaPj[y,{N', j'}, {N, j}, J] 


{n'{l> l/2)i'||AV^(9x)||nO l/2)j) 
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TABLE II: The same as Table [I]for the reduced matrix elements of the abnormal parity operators (see Eqs. (©). Recall that 
Q r f J (qx) and f2j J (gx) are not independent, see Eq. ([9}. 



Mathematica function 


Reduced matrix element 


DeltaJ[y, {n', j'}, {N, j}, j] 


(n'(l' l/2)j'||A* 3j (gx)||n(/ l/2)j) 


SigmaPj[y, {n\ j'}, {n, j}, j] 


(n'(l' l/2)j'||E: 7 M '(gx)||n(/ 1/2) j) 


SigmaPP j[y, {n', j'}, {N, j}, j] 


{n>(l> l/2)j'||E'; M ^gx)||n(/ 1/2) j) 


OmegaJ[y, {n', j'}, {N, j}, j] 


(n'(l' l/2)j'\\ny->(qx)\\n(l 1/2) j) 


OmegaP j[y, {n', j'}, {N, j}, J] 


(n'(l' l/2)f\\n'f 1J (^)\Hl 1/2) j) 



Notice that Mathematica has built in functions for the Wigner 3j and 6j symbols. For the practical reason of 
avoiding warning messages in the output, we explicitly imposed that a 3j symbol is put to zero when the conditions 
for its existence are not satisfied, while the built-in function is used otherwise. A similar procedure was used for 6j 
symbols. 

There is no built-in function for the Wigner 9j symbols in Mathematica, therefore we had to build one, using the 
definition: 

h I ju \ = y>D 2 *(2. 9 + 1) ( * * J ? } \ j * t \ { J i 3 J r i \ (46) 

Jl3 J24 J 



'I 



^ 1 J34 J 9 J 1 h 9 J24 J 1 9 ji j3 



The code implements an explicit check of the conditions to have non-zero reduced matrix elements: (i) parity 
conservation, (ii) that the nodal numbers n and n' are positive: n > 0, n' > 0, and (hi) that j,j',J satisfy the 
triangular inequality. Only if all three conditions are satisfied, the code proceeds with the evaluation of the function 
describing the reduced matrix element, otherwise the function is put to zero. 

In Appendix A we include the text of a Mathematica notebook (also available online [H) where SevenOperators 
used to calculate the first of the three examples we considered, 12 Q(v e , e~) 12 N(g.s.), yielding the result shown in Eq. 
©. 

We acknowledge support from the U.S. Department of Energy under grants DE-FG02-00ER-41132 and DE-FC02- 
01ER41187 (SciDAC), and from ORNL contract 4000058152. 
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Appendix A: Mathematica notebook calculation of 12 C(v e ,e ) 12 N(1 + 1) 

This notebook illustrates the uploading of the package (specifically, uploading the master file of the modular version) 
and the use of the seven basic operator matrix elements to construct the nuclear multipole operators in Sec. HVl The 
application is the charge-current cross section 12 C(v e , e~) 12 N(l + l), the first example discussed in this paper. The 
notebook describes points where a user will need to insert steps to account for the specific electroweak reaction, 
assumed one-body density matrix, isospin matrix elements, and allowed multipolarities. SevenOperators makes an 
otherwise tedious computation very immediate, eliminating much of the potential for algebraic errors. The final result 
in this notebook can be further manipulated to improve its form, but we have not done this here, in order to keep 
this appendix readable. Such final manipulations are also a matter of personal taste. 



Example of use of SevenOperators: calculation of neutrino charged current 
cross section on 12 C 
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. This notebook illustrates the use of the Mathematica package SevenOperators to calculate the 
neutrino charged current cross section for the reaction 12 C(v, e~) 12 N(g.s.), as given in Eq. (33) of this 
paper. The final result is Eq. (36). Although we do a specific example, the discussion here points out 
those steps where a user must supply information to modify this script for other possible transitions and 
reactions. 

First, we upload the package that computes the single-particle reduced space-spin matrix elements of 
the seven basic operators. The package must be placed in a directory searched by Mathematica, 
which one can determine with the command $Path. 



In[1]:= | << 7o_master_f ile.m 

. We input the onebody density matrix, 

which here is defined as the set of coefficients relating manybody matrix elements reduced in angular 
momentum and isospin to single nucleon matrix elements reduced in angular momentum and 
isospin. In this example calculation we truncate the density matrix expansion to the 1 p shell, 

so that there are four independent density matrix elements. We denote these as$ {1/1] , ${1,3}, 

${3,i } , and ${3,3} , and we use below the relation $ {3r 1} = ( $ + + $ ) / 2 and $ { 1r3] = ($ + - $ ) / 2, 

where $ + and $~ are the symmetric and antisymmetric combinations 
introduced in the paper because of operator time reversal properties. Because the 
transition being considered in this example picks out the J and T of the transition, 

it is sufficient to label the density matrix as $\2 S \2j\- Thus [${2j',2j), \N', j '), {N, /}} gives the density 
matrix element and the SevenOperators calling arguments. The $\2\\2 j\ enter symbolically , 

but could easily be replaced by numerical values. 



In[2]:= 



■ -{{*,i.i„{i, ^}}- 

{<#* ♦ r>/2,{i, i}, {i, \}}, 
{or - r>/2, {i, \), {i, \}}, 
{,„,,„ {l, |},{l, \}}}-, 



We have defined the density matrix in terms of doubly reduced 
to doubly reduced matrix elements ( which assumes good isospin labels), 

but our charged current cross section formulae {either Eq. (33) or the general case Eq. (24)) involve 
manybody matrix elements that are not reduced in isospin. Thus we can define anew density 
matrix that relates singly reduced manybody matrix elements to doubly reduced manybody matrix 
elements. The user in general will have to insert the proper Wigner Eckart 3 J symbol here, 

as it depends on the isospins of the initial and final states and the isospin of the transition 
(which could be isoscalar for electromagnetic or weak neutral current interactions). Our 
example has T f = M Tf = l, t = m t = l, and t± = m t . = o, 

so the Wigner Eckart factor {-if*-*** ThreeJSymbol[{T f , -M Tf ), {T, M T ], {7^, M T .}] = 1 / \[J , 

which we then fold into the density matrix. 



In[3]:= 



S T = {{b[[1, 1]]/ VT, B[[l, 2]], B[[l, 3]]}, 
[[2, 1]] /VT, B[[2, 2]], B[[2, 3]]}, 
[[3, 1]] / VT, B[[3, 2]], fi[[3, 3]]}, 
{ffi[[4, 1]] /VT, S[[4, 2]], iE[[4, 3]]}}; 



Using the package and the formulae in Eq. (30), 

we calculate the space - spin matrix elements of the operators M 5 , 

L5, T e/5 , T mag . As the density matrix in this example connects the 

manybody matrix elements to doubly reduced single nucleon matrix elements, 
the user has to multiply the package results for the seven basic operators by the appropriate 

single nucleon reduced isospin matrix element. Now for the charged current neutrino 

reaction under consideration, the isospin operator is t + = (r x +r 2 )l 2 = -r lrMT=+1 / ^[J , 
where r 1iMr=+1 denotes the M T = +1 component of the isospin spherical tensor operator 

(see the paper, if additional explanation is needed). We have already used the Wigner Eckart 
theorem to remove magnetic isospin quantum numbers. As < 1/2/1 r z // 1/2 > =^J~6 , 
we see that the single nucleon isospin factor is - y/~3, which we insert by hand 
below. Users will need to supply analogous factors, depending on the process 
(charged current weak, neutral current weak, electromagnetic) being considered. 



ln[4]:= M 5 [y_, {Np_, jp_}, {N_, J_] 
I*q 



ln[5]:= 



^5 |f a [1] OmegaPJ[y, {Np, jp} , {N, j}, J] + 1 wF P [l] SigmaPPJ[y, {Np, jp}, {N, j}, J] j (-VT) 
L 5 [y_, {Np_, jp_}, {N_, j_}, J_] := I |f a [1] -I_2.F P [l]j SigmaPPJ[y, {Np, jp}, {N, j}, J] (-VT) 




. Now, we fold the density matrix with the corresponding doubly-reduced single nucleon matrix elements. 



In[8]:= 



ln[9]:= | SummedL 5 [y_] : = Sum[S T [ [i, 1]] *L 5 [y, J! T [[i, 2]], 9> T [[i, 3]], 1], {i, 1, Length [ ffi T ]} ] 
ln[10]:= I SutranedT al5 [y_] : = Sum[$ T [ [i, 1] ] * T el5 [y, & T [ [i, 2] ] , 5 T [ [i, 3] ] , 1] , {i, 1, Leng 



ln[11]:= | SuiranedT mag [y_] : = Sum[* T [ [i, 1] ] * T mag [y, 2]], B T [[i, 3]], 1], {i, 1 , Length [ffi T ] } ] 

. We calculate the three terms governing the nuclear response (see Eq. (33)), removing the common 
term e~ 2y . 



In[12]:= 



ln[13]:= 



FirstTerm[y_] := Simplif y ^Exp [ 2 y ] fsummedM 5 [y] + — SummedL 5 [y] J Conjugate ^SummedM 5 [y] + — SummedL 5 [y] j 



{M, q, y, 01, F P [1] , Fi[l] , F A [1] , , lit* , , * ( 3,3}} e Reals 



SecondTerm[y_] := 

Simplify [Exp [2 y] ( SummedT e i 5 [y] Conjugate [SummedT e i 5 [y] ] + SummedTmag [y] Conjugate [SummedT mag [y] ] ) 




ln[14]:= I ThirdTerm[y_] := Simplif y [Exp [2 y] Re [ SummedT mag [y] Conjugate [SummedT e i 5 [y] ]] , 




Now, one can compose the full expression of the cross section for 12 C(v, e ) 12 N(g.s.), reinserting e 2y . 



In[15]:= 



2e"2 

Xsection[s_, y_] :- G F * 2 Cos[e c ] '2 F el Cos[S/ 2] "2 Exp[-2 y] 



In A 2 

FirstTerm[y] + I + Tan[e/2]~2| SecondTerm [y ] 

2q"2 




2 Tan [6/ 2] 



+ Tan[S/2]"2 ThirdTerm[y] 
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The cross section is printed. 



In[16]:= 



Out[16]= 



-2e"' y E'Cos - Cos [9 ] F el G 2 

71 L 2 J 



— q (5 (-l + 2y) +V2" (VF (5-4y) 1^3,3, + 5 (-2 + y) f) ) Tanf-1 -— + Tan[-| 

25 M 71 * 1 2 V q 2 

(5 *>(!,!) ((-l + 2y)M[l] + 2 Fi[l]) +VT (5r ((-2+y)(/[l] +Fl[l]) + V5"<|r {3 ,3, ( (5 - 4 y) (7 [ 1 ] + 5 F, [ 1 ] ) ) 

1 I il re, 1 

Fa [ 1 ] + + Tan - 

900 71 [ 2 q 2 L 2 J 

V2" (5*- ( (-2 +y) u[l] + Fi[l]) + V5> { 3, 3 ) ((5-4y) <i[l] + 5F![1]))) 5 
4 (5 (-l + 2y) |f-„, 1)+ Vr (VF (5-4y) ^3,3) + 5 (-2+y) r)) 2 F A [l] 2 ) + 



-q 2 5Vu,i) ((-l + 2y) f7[l] +2Fi[l]) 



900 M 7Tq 2 

(5 -JY (l + 2y) u * {1>1) (HF,[1] -q 2 F p [l]) + 2 VF (-5 + 2 y) w i/r {3 , 3) (MF A [1] - q 2 F p [ 1 ] ) 



5 (3 q 2 r F A [1] + 4 (-1 +y) i)f (MF A [1] - q 2 F P [1] ) ) ' 



• One can continue to manipulate this expression, according to one's taste. It is straightforward to show 
this result is equivalent to Eq. (36). 



